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Dissipationless collapses in Modified Newtonian Dynamics (MOND) have 
been studied by using our MOND particle-mesh N-body code, finding that the 
projected density profiles of the final virialized systems are well described by 
Sersic profiles with index m<4 (down to m ~ 2 for a deep-MOND collapse). 
O ' The simulations provided also strong evidence that phase mixing is much less 

effective in MOND than in Newtonian gravity. Here we describe "ad hoc" 
numerical simulations with the force angular components frozen to zero, thus 
producing radial collapses. Our previous findings are confirmed, indicating that 
possible differences in radial orbit instability under Newtonian and MOND 
1 gravity are not relevant in the present context. 
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1. Introduction 

In the Lagrangian formulation of Milgrom's Modified Newtonian Dynamics 
(MOND)^'^ the Poisson equation V^(/)n = inGp for the Newtonian poten- 
tial (/)N is replaced by the field equation for the MOND potential (f> 

V-[/.(||V0||/ao)V0] =47rGp, (1) 

where p is the density distribution, gq 2± 1.2 x 10~"'^°ms~^ is a character- 
istic acceleration, ||...|| is the standard Euclidean norm, and in finite mass 
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systems V0 for ||x|| oo. The MOND gravitational field experienced 
by a test particle is g = — V0, and ^ y for y <C 1 and ^ 1 for y S> 1 
(typically /i = y / ^/T+lp). In the 'deep MOND regime' describing low- 
acceleration systems (||V0|| ^ ao, hereafter dMOND), = y and so 
Eq. (1) simplifies to V • (||V^||V^) = AwGaop- The source term in Eq. (1) 
can be eliminated by using the Poisson equation, giving 

/i(||V0||/ao)V0 = V0N + S, (2) 

where S is a solenoidal field dependent on p and in general different from 
zero; when S = Eq. (2) can be solved explicitly in terms of V(/)n- This 
reduction would be most useful for numerical simulations. Unfortunately 
S = only for very special^'* configurations. In addition, though the field 
S has been shown to be small^'^ for some configurations, neglecting it when 
simulating time-dependent dynamical processes has dramatic effects such 
as non-conservation® of the total linear momentum. 

Several astronomical observational data appear consistc^nt^'^ with the 
MOND hypothesis, and also a relativistic version'^ of MOND is now avail- 
able, making it an interesting alternative to the cold dark matter paradigm. 
However, dynamical processes in MOND have been investigated very lit- 
tle^°~^^ so far, mainly due to diflaculties posed by the non-linearity of 
Eq. (1). In a recent paper^ (hereafter NLC) we presented the results of 
N-body simulations of dissipationless collapse in MOND obtained with our 
N-body code which solves Eq. (1) exactly. In particular, we obtained clear 
indications that phase mixing is much less effective in MOND than in New- 
tonian gravity. Here, after summarizing the main results of NLC, and re- 
stricting for simplicity to the Newtonian and dMOND regimes only, we 
address the problem of the importance of the force angular components in 
the relaxation process, by running "ad hoc" numerical simulations with the 
force angular components frozen to zero. 

2. The N-body code and the numerical simulations 

Our MOND N-body code^'^ is based on a particle-mesh scheme with 
quadratic spline interpolations. The spherical grid on which Eq. (1) is solved 
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is made of Nr x x iV^ points. We use leap-frog time integration, where 
the adaptive time-step is the same for aU particles. All the computations 
on the particles and the particle-mesh interpolations can be split among 
different processors, while the iterative potential solver computations are 
not performed in parallel: however, at each time step we can use the po- 
tential previously determined as seed solution. We succesfully tested the 
code by running Newtonian simulations (i.e., by solving Eq. [1] in the limit 
/X = 1), and comparing the results with those of simulations performed with 
our FVFPS^^'^® treecode starting from the same initial conditions. We also 
verified that the code reproduces the Newtonian and MOND conservation 
laws. In fact, 2K + W = for virialized systems in MOND and in New- 
tonian gravity, where K is the total kinetic energy and W = Tr Wij is the 
trace of the potential energy tensor^''' 



Note that in MOND K + W is not the total energy, and is not con- 
served. However, W is conserved in the limit of dMOND, being W{t) = 
-(2/3)v/GaoM3 for all systemsi'i^-^i ^^^^-^ ^0^^! j^^ss M,. 

The choice of appropriate scaling physical units is an important as- 
pect of MOND N-body simulations. A full discussion of this point can be 
found in NLC; here we just recall that, while in Newtonian simulations 
the natural scaling units are = r*^^(GM,)~-^/^, = (GM*)^/^r7^^^, 
and E^n = GM^'^r''^, in the dMOND case one has t^^ = r*(GM*ao)~^'''^, 
v^d = (GM^ao)^/^, and E^^ = {GaoY^^M^^^"^ , where and are the 
length and mass units in which the initial conditions are expressed. In NLC 
we performed a set of N-body simulations of dissipationless collapses, start- 
ing from the same phase-space configuration, i.e., from a cold {2K/\W\ = 0) 
Plummer^^ sphere of total mass M* and "core" radius r*. The gravitational 
potential is Newtonian in simulation N and dMOND in simulation D; the 
full MOND simulations are not described here (see NLC). All the simula- 
tions (iV = 10^ particles, A'^,. = 64, N^} = 16 and A^,^ = 32) are evolved up 
to t = 150tdyn, where tdyn is defined as the time at which 2K/\W\ reaches 
its maximum value (tdyn ~ 2t*d in simulation D, and ~ 2t*n in N). The 
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Fig. 1. Time evolution oi2K/\W\, W, and K + for a dMOND (D) and a Newtonian 
(N) simulation. K, W, and K + W axe in units of E*d (D), and E^n (N). Note the 
long-lasting oscillations of the virial ratio in D. 



center of mass, and the modulus of the total angular momentum (in units 
of r«M*f*n for model N and of r*M*z;*d for model D) oscillates around 
zero with r.m.s <0.1r* and <0.02, respectively. The quantities K + W in 
simulation N, and W in simulation D are conserved to within 2% and 0.6%, 
respectively (see Fig.l). The final angle-averaged density profiles are fitted 
with a 7-modeP^'^^ 

(3 - 7)M,re 



pir) = 



(4) 



for each end-product we also measure the axis ratios c/a and b/a of the 
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inertia cllipsoid^^'^^ {a >b> c), the ellipticity e of the projections along the 
principal axis, and the corresponding circularized effective radius R^. The 
resulting circularized projected density profiles are fitted with the Sersic^'''^* 
law 

/(i?) =/(i?e)e-''('")[(«/«")''"""i], (5) 

where^^ b{m) ~ 2m — f/3 + 4/405m. Note that m is the only free parameter, 
because i?e and /(-Re) are fixed by particle count. 

3. Results 

3.1. Structure and kinematics of the collapse end-products 

All the final virialized systems in NLC depart significantly from spherical 
symmetry. In particular, the D end-product is triaxial (c/a ~ 0.2, b/a 0.4; 
0.5 < e < 0.8), while model N is oblate (c/a ~ c/b ~ 0.5; < e < 0.5). These 
values are consistent with those observed in real ellipticals, with the excep- 
tion of model D, which would correspond to an E8 galaxy. Thus, MOND 
gravity could be able to produce some system that would be unstable in 
Newtonian gravity. The final angle-avcragcd density of model N is well de- 
scribed by Eq. (4) with 7 ~ 1.7, while 7 ~ in model D. Equation (5) 
fits remarkably well the final surface density profiles (best-fit index m ^ 4 
in model N, and m 2 in D, see Table 2 in NLC), with average residu- 
als 0.05 <(AS'B)< 0.2, where SB = -2.51og[/(i?)/7(iie)]. Note that the 
fitting radial range 0.1 ^ R/Rg ^ 10 is comparable with or larger than 
the typical ranges^° spanned by observations. The internal kinematics of 
the end-products is quantified by the angle-averaged radial and tangen- 
tial components of their velocity-dispersion tensor (cr^ and at), and by the 
anisotropy parameter (3{r) = 1 — Q.5a^ / . All systems are strongly radially 
anisotropic outside the half- mass radius. The ar profile decreases steeply in 
the final state of model N, while it presents a hole in the inner regions of 
model D. In addition, model D is radially anisotropic {(3 0.4) even in the 
central regions, where model N is approximately isotropic {f3 ~ 0.1). The 
line-of-sight velocity dispersion crios declines steeply within R^ in model N, 
while the D profile is significantly flatter. 
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3.2. Phase-space properties 

In Newtonian gravity, collisionless systems virialize through violent relax- 
ation in few dynamical times, as predicted by the theory^^ and confirmed 
by numericaP'*''^^ simulations. Due to the non linearity of the theory, the 
details of relaxation processes and virialization in MOND are much less 
known. In Fig. 1 we show the time evolution of 2K/\W\, W, and K + W 
for simulations D and N of NLC. In simulation N, 2K/\W\ has a peak, then 
oscillates, and eventually converges to the equilibrium value 1; the total en- 
ergy K + W is nicely conserved. The time evolution is significantly different 
in simulation D, where W is constant as expected, but 2K/\W\ still oscil- 
lates at very late times because of the oscillations of K. A different view of 
phase space is given in Fig. 2, whore we plot time snapshots of the particles 
radial velocity vs. radius for simulations D and N. At t = tdyn (time of the 
peak of 2i4r/|M^|), sharp shells in phase space are present, indicating that 
particles are moving in and out collectively and phase mixing has not taken 
place yet. At significantly late times {t = 44tdyn), when the systems are 
almost virialized (2i^/|W| ~ 1), phase mixing is complete in simulation N, 
but phase space shells still survive in model D. Finally, the bottom panels 
show the (r,Vr) plane at equilibrium {t = ISOtdyn), when phase mixing is 
completed also in model D: note how the populated regions are significantly 
different for the two models. Thus, NLC results indicate that phase mixing 
in more effective in Newtonian gravity than in MOND. Here we address 
the issue of the importance of the force angular components during the col- 
lapse. In fact, one could speculate that a different behavior of radial orbit 
instability in MOND and in Newtonian gravity could be at the origin of 
the different time scale of phase mixing. As shown in Fig. 3, also the new 
simulations confirm that phase mixing is less effective in MOND than in 
Newtonian gravity (even though virialization times are now longer than in 
NLC models, due to the reduced number of active degrees of freedom). 

In NLC we obtained additional information on the relaxation process 
from the differential energy distribution^'' N{E) (i.e. the number of particles 
with energy per unit mass between E and E+dE). In Newtonian gravity (j) is 
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Fig. 2. Radial-velocity vs. radius of 32000 particles randomly extracted from simula- 
tions D and N. Vr is in units of ii«d (D) or u«n (N). 



usually set to zero at infinity for finite-mass systems, so E = u^/2-|-(^(x) < 
for bound particles; in MOND all particles are bound independently of their 
velocity, because (j) is confining, and all energies are admissible (see Fig. 5 
of NLC). Given that the particles are at rest at t = 0, the initial N{E) de- 
pends only <^(x) at t = 0, and it is significantly different in the Newtonian 
and MOND cases. In accordance with previous studies, we found that in 
the Newtonian case the final N{E) is well represented by an exponential 
function^^'^^'"'^"*^^ over most of the populated energy range. In contrast, in 
model D the final N[E) decreases for increasing energy, qualitatively pre- 
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serving its initial shape. We interpret this result as another manifestation 
of a less effective phase space reorganization in MOND than in Newtonian 
collapses. 
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Fig. 3. Phase space sections for simulations with frozen angular force components (bot- 
tom panels), and with the three force components active (top panels). 



4. Conclusions 

We presented results of dissipationless collapses in MOND, focusing on the 
relaxation process. The main results can be summarized as follows: 

1) Newtonian collapses produced cuspier density profiles than MOND 
simulations (7 ~ 1.7 and m 4 vs. 7^0 and m ~ 2). In both cases 
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the Sersic fits are remarkably good over a large radial interval. In addi- 
tion, Newtonian models arc isotropic in the central regions, while dMOND 
models arc radially anisotropic down to the center. 

2) In NLC we found that final states of full MOND models, if interpreted 
in the context of Newtonian gravity, are characterized by a dividing radius 
of the order of -Re, separating a baryon-dominated inner region from a 
dark-matter dominated outer region, in accordance^^'^'' with observations 
of elliptical galaxies. However, we were not able to reproduce the observed 
scaling laws of elliptical galaxies under the assumption of a luminosity- 
independent stellar mass-to-light ratio. 

3) Phase mixing is less effective (and stellar systems take longer to 
relax;) in MOND than in Newtonian gravity. This behavior is confirmed by 
numerical simulations in which the angular force components are frozen to 
zero, so possible differences in radial orbit instability between Newtonian 
and MOND gravity are not relevant in the present context. Our results on 
mixing suggest that merging could take longer in MOND than in Newtonian 
gravity; on the other hand, analytical estimates of the two-body relaxation 
time seem to indicate"^* the opposite, predicting shorter dynamical friction 
time-scales in MOND than in Newtonian gravity: the next application of 
our code will be the study of galaxy merging in MOND. 
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